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ABSTRACT 


An  algorithm  is  developed  which  enables  the  computation 
of  horizontal  range  and/or  depth  from  a  submerged  sound 
source,  using  ray  acoustics  and  the  Lloyd's  mirror 
interference  effect.  The  solution  is  based  on  Snell's  law 


and  involves  integrating  multipath  sound  rays  to  find  the 
difference  in  length  between  the  direct  and  surface 
reflected  sound  paths  from  the  source  to  the  receiver.  This 
difference  in  path  length  is  directly  related  to  the 
observed  Lloyd's  mirror  interference  pattern. 

No  assumptions  as  to  the  mathematical  characteristics  of 
the  sound  velocity  profile  (SVP)  are  made  nor  are  far  field 
approximations  used.  The  solution  is  as  accurate  as  the  SVP 
data  input  to  the  problem. 

A  general  computer  flow  chart  and  basic  language  program 


are  provided  to  allow  local  commands  the  capability  of  pre¬ 
mission  planning  based  on  specific  operating  area 
environmental  information. 
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I.  INTRODUCTION 


Sound  propagation  between  a  source  and  a  receiver  in  the 
ocean  may  occur  by  many  different  ray  paths,  each  dependent 
on  existing  physical  and  acoustical  conditions.  Crucial  to 
the  determination  of  which  sound  paths  cure  of  importance  in 
any  given  problem  are  factors  such  as  the  sound  velocity 
profile,  the  depth  of  the  water,  the  bottom  contour  and 
composition,  the  sea  surface  conditions  and  the  intensity 
and  frequency  of  the  sound. 

This  thesis  deals  with  two  specific  paths  of  sound 
propagation  from  a  source  to  a  receiver.  The  first,  called 
the  direct  path,  does  not  touch  the  surface  or  bottom  of  the 
sea  as  it  travels  from  the  source  to  the  receiver.  The 
second  path  reflects  off  the  sea  surface  before  reaching  the 
receiver  and  is  called  the  reflected  path. 

The  following  explanation  of  the  reception  of  sound  from 
a  single  source  via  these  two  distinct  paths  is  taken 
directly  from  Urick:  [Ref.  1:  pp.  120-123] 

"If  the  sea  surface  were  perfectly  smooth,  it  would 
form  an  almost  perfect  reflector  of  sound.  The  intensity 
of  sound  reflected  from  the  smooth  sea  surface  would  be 
very  nearly  equal  to  that  incident  upon  it." 

He  continues: 

"A  criterion  for  the  roughness  or  smoothness  of  a 
surface  is  given  by  the  Rayleigh  parameter,  defined  as 
R*  kH  sinQ  ;where  k  is  the  wave  number,  (2tt)/a  (  A  is  the 
wave  length  of  the  sound  incident  on  the  surface),  H  is 


the  rms  (sea)  wave  height  (crest  to  trough),  and  &  is  the 
grazing  angle". 


Figure  1-1.  Rayleigh  Parameter  (R)  for  Describing  Sea 
Surface  Roughness. 

"When  R  <  <  1 ,  the  surface  is  primarily  a  reflector  and 
produces  a  coherent  reflection  at  the  specular  angle  equal 
to  the  angle  of  incidence.  When  R>>1,  the  surface  acts 
as  a  acatterer,  sending  incoherent  energy  in  all 
directions.  When  certain  theoretical  assumptions  are 
made,  (See:  P.  Beckmann,  and  A.  Spizzichino:  Scattering  of 
Elec tromagnetic  Waves  From  Rough  Surfaces,  p.  93,  The 
Macmillan  Company,  N.Y.,  1963.)  the  (pressure)  amplitude 
reflection  coefficient  m  of  an  irregular  surface  (defined 
as  the  ratio  of  the  reflected  or  coherent  pressure 
amplitude  of  the  return  to  the  incoherent  pressure 
amplitude)  can  be  shown  to  be  given  by  the  simple 
expression  m  «  exp  (-R)." 

Urick  goes  on  to  explain  the  image  interference. 

"When  the  sea  surface  is  not  too  rough,  it  creates  an 
interference  pattern  in  the  underwater  sound  field.  This 
pattern  is  caused  by  constructive  and  destructive 
interference  between  the  direct  and  surface  reflected 
sound  and  is  called  the  Lloyd  mirror,  or  image 
interference,  effect." 

In  Crick's  sound  propagation  in  the  sea  [Ref.  2:  p.  9-5],  he 

continues  the  discussion  of  surface  reflection  interference: 

"This  sound  field  may  be  divided  into  three  parts,  1)  a 
NEAR  FIELD  close  to  the  source  in  which  the  image  source 
is  too  far  away,  and  the  reflected  sound  is  too  weak  to 
produce  appreciable  interference,  2)  an  INTERFERENCE  FIELD 
in  which  there  are  strong  loops  and  nulls  in  the  signal 
received  by  a  receiver  moving  outward  in  range,  and  3)  a 
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FAR  FIELD  in  which  there  is  an  increasingly  out -of -phase 
condition  between  source  and  image  and  the  intensity  falls 
off  as  the  inverse  fourth  power  of  the  range  (transmission 
loss  increases  as  40  log  r).” 

The  interference  effects  felt  at  the  receiver  is 
expressable  as  the  ratio  of  the  combined  direct  and 
reflected  sound  intensities  to  the  direct  sound  intensity 
alone. 

Stated  mathematically: 


Intensity  Ratio 


1 1  ot  al 


I(direct+  reflected) 


or  more  precisely 
Intensity  Ratio  = 


Idirect 


Idirect 


2m 


rd\ 

. .  1  cos(wr) 

rr/ 


Where  RD  is  the  direct  ray  path  length,  RR  is  the  reflected 
path  length,  m  Is  the  surface  pressure  reflection 
coefficient,  u>  (where  «  =  2xf)  is  the  angular  frequency  of 
the  sound,  and  r  is  the  time  delay  caused  by  the  difference 
in  the  two  path  lengths.  A  formal  derivation  of  the 
intensity  ratio  is  given  in  Appendix  A.  A  more 
comprehensive  treatment  of  the  effects  of  surface  roughness 
on  acoustic  waves  is  given  by  Clay  and  Medwin  (Ref.  3]. 
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II.  APPLICATION  of  LLOYD  MIRROR  INTERFERENCE  PHENOMENON  to 
LOCATE  the  HORIZONTAL  RANGE  and  DEPTH  of  the  SOUND 


SOURCE  FROM  the  RECEIVER 


Tactical  tracking  information  from  the  Lloyd's  mirror 
phenomenon  has  been  used  operationally  by  the  Navy  for  many 
years.  In  an  excellent  report  on  the  analysis  of  Lloyd's 
mirror  phenomenon  under  constant  velocity  gradient 
conditions,  Shudde  [Ref.  4:  p.  3],  discusses  the  difficulty 
of  solving  the  surface  reflection  interference  problem. 
Accurate  error  estimation  for  any  acoustical  condition  other 
than  the  isovelocity  case  is  not  practical  unless  a  large 
computer  is  used.  Since  Shudde's  1977  report,  the 
advancements  in  memory  size,  and  speed  of  operation  of  small 
portable  desk  top  computers  now  allow  an  accurate  solution 
of  source  position  for  all  sound  velocity  conditions.  The 
following  mathematical  solution  is  based  on  the  lengths  of 
the  sound  paths  determined  from  Snell's  Law.  As  explained 
by  Urick,  [Ref.  1:  p.  116], 

"One  of  the  most  important  practical  results  of  ray 
theory  is  Snell's  Law,  which  describes  the  refraction  of 
sound  rays  in  a  medium  of  variable  velocity." 

The  solution  is  not  dependent  on  any  isovelocity  or  constant 

gradient  conditions  other  than  the  reasonable  assumption 

that  no  horizontal  sound  speed  gradients  exist  throughout 

the  operating  area.  The  mathematical  solution  presented  in 

Section  III  makes  no  assumption  as  to  the  mathematical 


characteristics  of  the  sound  velocity  profile  nor  does  it 
use  far  field  approximations.  Past  solutions  depended  on 
specific  sound  velocity  profiles,  such  as  an  isovelocity  or 
"straight  line"  profile  which  resulted  in  significant 
positional  errors  for  some  tactical  applications  when  the 
actual  sound  velocity  profile  did  not  match  the  predicted 
model. 

In  applying  the  solution  to  the  Lloyd's  mirror 
interference  effect  one  must  understand  the  relationship  of 
the  propagating  sound  to  the  interference  pattern.  Coppens, 
Sanders  and  Dahl  [Ref.  5  :  p.  55  ]  explained, 

"...the  noise  radiated  by  underwater  vehicles  consists 
of  two  contributions:  broadband  noi se  (similar  to  the 
hissing  found  between  FM  radio  stations)  and  tonals 
(single  frequency  tones,  like  whistles).  Thus,  the 
radiated  noise  consists  of  all  frequencies  between  certain 
broad  limits  with  certain  specific  frequencies  occurring 
with  high  intensity.  Each  of  these  frequencies,  whether 
of  broadband  noise  or  tonal,  has  its  own  interference 
pattern. 

Because  each  frequency  has  its  own  period  and 
wavelength,  each  will  have  a  different  phase  delay  for 
propagation  over  the  direct  and  surface  reflected  paths. 
As  a  result,  the  ranges  between  source  and  receiver  for 
which  the  surface  reflected  signal  tends  to  cancel  the 
direct  signal  will  be  different  for  each  different 
frequency.  Thus,  the  ranges  at  which  different 
frequencies  will  have  nulls  at  the  receiver  will  be 
different.  This  means  that  as  the  range  between  source 
and  receiver  changes,  the  frequencies  which  have  nulls  at 
the  receiver  will  also  change." 

Figure  II-1  shows  a  schematic  time  vs.  frequency  plot  of  the 

interference  pattern  as  the  distance  between  the  source  and 

receiver  is  decreased  and  then  increased.  The  distinct 

tonals  will  most  readily  show  the  effects  of  Doppler  shift 


(apparent  frequency  change)  as  minimum  separation  range,  is 
passed  (see  t  >  0  in  Figure  II-1).  The  light  (unshaded) 
hyperbolic  shaped  regions  centered  about  t  «  0  are  known  as 
the  Lloyd's  mirror  interference  pattern  for  the  broad  band 
noise.  The  frequency  values  Ff ,  F2 ,  F3 ,  F4  represent  the 
minimum  frequency  of  each  hyperbola.  The  actual  frequency 
from  the  sound  source  that  is  being  nullified  by  the  surface 
reflection  interface  varies  with  time  as  the  hyperbolas  are 
traced  out.  The  measured  separation  in  hertz  between  any 
two  successive  nulls  at  any  time  t,  is  termed  the  difference 
frequency,  A  f ,  and  is  the  frequency  value  required  in  the 
problem  solution  discussed  in  the  following  Section.  (See 
Section  I1I.D.2) 


III.  RANGE  VS.  DEPTH  SOLUTION  USING  the  LLOYD'S  MIRROR 

INTERFERENCE  PHENOMENON 


A.  THEORETICAL  BACKGROUND 
1 .  Snell's  Law 

Snell's  Law,  [Ref.  1:  p.  116,  Ref.  6:  p.  120],  is  a 
result  of  Fermat's  minimum  time  principle  which  states  that 
each  ray  will  travel  the  path  requiring  minimum  time  of 
travel  from  source  to  receiver. 

Mathematically,  Snell's  Law  is  expressed  as: 

sin  0 

Snell  Constant  *  - — - - 

c(z) 


Where  0  is  the  angle  that  the  ray  makes  with  the  vertical 
and  c(z)  is  the  sound  velocity  for  a  given  depth  z.  (See 
Figure  III-1 ) 


Rinitial 


Surface 


zinitial' 


origin 


Source 


+  Range, 
Meters 


Raypath  direction  at  Depth  Z±, 

Range  R^ 


Depth  ( z ) , 
Meters  + 


Figure  III-1.  Snell's  Law  Applied  to  Initial  Sound  Source 

Ray  Angle. 
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The  sound  ray  maintains  the  same  Snell  Constant  value 
throughout  its  transmission  path.  Sound  velocity  changes 
result  in  an  increasing  or  decreasing  angle  0. 

2.  Rav  Path  Integration 

Expanding  the  diagram  of  the  sound  ray,  the 
differential  raypath  length,  ds,  can  be  described  as  a 
function  of  the  differential  change  in  depth,  dz,  and  the 
differential  change  in  horizontal  range,  dr.  (Figure  III-2) 


Sound 

Ray 

* 

✓ 


ds  = 


♦  dr2 


Figure  III-2.  Pythagorean  Theorem  Applied  to  Ray  Path 

Length. 


To  calculate  the  actual  path  length  of  the  sound 
ray,  smaller  and  smaller  segments  of  path  length,  ds, 
horizontal  range,  dr,  and  depth,  dz  are  used  until  the 
differential  segments  ds,  dr,  and  dz  are  approached.  This 
tends  to  smooth  the  curve  of  the  ray  into  a  sequence  of 
straight  line  segments.  The  path  length  traveled  by  a 
particular  ray  can  then  be  determined  by  integrating  between 
the  source  depth  and  receiver  depth.  The  possibility  of  the 
direct  ray  crossing  the  same  depth  more  than  once  is  not 
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considered  a  problem  in  this  solution  due  to  the  short 
ranges  involved  in  Lloyd's  mirror  observations  and  the  long 
ranges  required  for  a  sound  ray  to  curve  back  through  the 
same  depth  before  reflecting  off  the  surface.  To  handle  the 
fact  that  the  surface  reflected  ray  does  cross  the  same 
depth  twice  in  the  ranges  pertinent  to  this  discussion  and 
to  avoid  errors  in  integration,  the  concept  of  a  geometric 
equivalent  image  [Ref.  1:  p.  123],  of  the  source  is  used. 
The  image  reflected  path  and  the  actual  reflected  path  are 
the  same  length.  (Figure  III-3) 

/b'i  Image  Reflected  Path 

"2b  K 


c(z)  \  Surface  o 
SVP  — 1 


+  Range 
Meters 


Source 

♦*b' 


Actual 

Reflected 

Path 


+ Depth l 
Meters 


Figure  III-3.  Actual  and  Image  Sound  Paths  From  Source  to 
Receiver  for  a  Given  Sound  Velocity  Profile 
(SVP). 


The  curved  sound  path  shown  in  Figures  III -3  and 
III-4  is  due  to  the  change  (in  general)  in  the  sound 
velocity  as  the  depth  changes.  It  is  possible  to  show  that 
the  path  from  A  to  B  is  the  same  length  as  the  path  from  A 


I 

to  B  .  Also,  by  Inspection  of  the  sketch,  the  horizontal 

i 

rang*  from  A  to  B  or  A  to  B  is  the  same.  Therefore,  for 
the  reflected  path  calculations  of  path  length  and 
horizontal  range,  the  image  reflection  depth  of  -z^  will  be 
used. 

B.  HORIZONTAL  RANGE  SOLUTION 
1 .  Direct  Path 

The  horizontal  range  traversed  by  the  direct  path 
sound  ray  can  be  calculated  from  Snell's  Law  and  the 
Pythagorean  Theorem  for  right  triangles  applied  to  Figure 
III-2. 

From  (?)  ,  the  Snell  Constant,  K-j  is: 
sin  0* 

*1  .  - 

cCzj^) 


This  is  the  Snell  constant  for  the  particular  direct  path 
ray  that  ultimately  travels  from  the  source  to  the  receiver. 


but  initially  starts  at  depth  and  moves  in  direction  6^. 
Using  the  Pythagorean  Theorem  on  Figure  III-2 


sin  0 


Solving 


dr 


(dz)2 


(1A)  for  sin  0, 


sin  0  ■  K-j  c(z) 


where  0  »  variable  for  the 
particular  direct  path  ray 
tagged  with  Snell  constant 

K 1 » 
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sotting  (2)  and  (3)  equal  to  each  other  yields: 


dr 


The  right  side  of  @  can  be  simplified  by  dividing 
numerator  and  denominator  by  dr. 


Solving  for  dz/dr  follows  in  steps: 


20 


Now,  solving  for  dr,  the  incremental  horizontal 


range  yields: 


dr 


Integrating  @  with  respect  to  depth,  with  the 
limits  of  integration  being  the  depth  of  the  source  and  the 
depth  of  the  receiver,  yields: 


R  (direct  path) 


Source 

Depth 


Actual 

Receiver 

Depth 

(along  path  of  ray 
"tagged"  K-j ) 


The  "Range"  is  always  positive  and  thus  the 
appropriate  +  or  -  sign  is  used  to  make  this  integral  always 
positive  regardless  of  the  integration  limits.  For  ease  of 
discussion  in  the  rest  of  this  solution,  the  receiver  will 
be  assumed  to  be  at  a  shallower  depth. 
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2.  Reflected  Path 

Recalling  the  inage  depth  discussion  and  noting  a 
different  angle,  ©ir  for  the  reflected  ray  is  required 
(Figure  III-4),  which  results  in  a  different  Snell  constant, 
K2,  where 

sin  (8j) 

K2  »  - 

c(z± ) 

It  follows  directly  that  (7)  can  be  rewritten  for  the 
reflected  horizontal  range  as: 


R  (reflected  path) 


Source 


Image 

Receiver 

Depth 
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Which  can  be  simplified  by  multiplying  by 


as  follows: 


Receiver 

Depth 


2 


dz 


Substituting  the  equivalent  of 


from  the  reciprocal  of 


Source 

Depth 


S  (Direct  Path) 


1  + 


Receiver 

Depth 


K.,2  c  ( z ) 2 
1  -  2  c  ( z ) 2 


dz 


Source 

Depth 


S  (Direct  Path)  ■ 


Receiver 


c(z)2  +  Ki 2 
c  ( z )  ]  2 


c(z)2 


dz 


!  Source 

Depth 

© 


Actual 

Receiver 

Depth 


S  (Direct  Path) 


dz 


/ 


1  -  [K-,  c ( z )  ] ‘ 


2.  Reflected  Path 

Using  the  same  arguments  about  the  reflected  path  as 
presented  for  horizontal  range  and  substituting  K2  for 
and  image  receiver  depth  into  (Tl)  yields: 


Source 

Depth 


S  (Reflected  Path) 


1 


\Jl  -  [K2  c ( z )  ] 2 


dz 


Image 

Receiver 

Depth 


D.  EXPLANATION  OF  SOLUTION 
1 .  Some  Observations 

Equations  (f)  and  (12)  are  path  integrals  for  the 
direct  ray,  using  the  same  Snell  constant  and  integration 
limits.  Further,  equations  (?)  and  (?3)  are  path  integrals 
of  the  reflected  ray  and  thus  use  a  different  Snell  constant 
from  the  direct  path  integrals,  but  use  the  same  integration 
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limits  respectively  for  the  horizontal  range  and  ray  path 
lengths . 

2.  Solution 

The  key  to  the  eventual  numerical  solution  is  the 
equality  between  the  direct  path  horizontal  range  and  the 
reflected  path  horizontal  range.  Although  the  actual  path 
length  of  the  reflected  ray  is  longer  than  the  direct  path 
ray  between  common  source  and  receiver,  the  horizontal 
distance  traveled  by  both  is  the  same. 

By  calculating  all  of  the  initial  ray  angles  that 
are  possible  for  each  depth  of  the  source  of  interest,  two 
horizontal  range  arrays,  one  for  the  direct  path  and  one  for 
the  reflected  path,  are  formed.  Each  of  these  arrays  is 
linked  respectively  to  a  direct  and  reflected  ray  path 
length  array  by  the  corresponding  Snell  Constant  used  to 
determine  both  the  horizontal  range  and  sound  ray  length  for 
each  path,  (See  Figure  III-4).  For  example,  suppose  the 
initial  direct  path  ray  angle  of  15  degrees  from  the 
vertical  covered  a  horizontal  range  of  20  meters  in 
traveling  from  the  source  to  the  receiver.  Similarly,  the 
initial  reflected  path  ray  angle  of  8  degrees  from  the 
vertical  covered  the  same  horizontal  range  of  20  meters  from 
the  source  to  the  receiver.  Then  the  difference  between  the 
reflected  ray  path  length  for  8  degrees  and  the  direct  ray 
path  length  for  15  degrees  could  be  computed  (Figure  III-5). 


ii 


20  meters 


Range 

Meters 


Receivers 


8  degrees 
15  degrees 


.Source 


+  Depth 
Meters 


Figure  III-5.  Example  of  Reflected  Path  Ray  Interfering 
With  Direct  Path  Ray  at  the  Receiver. 


It  is  this  difference  in  ray  path  length, 
AS  «  -  Sp,  for  a  wide  band  of  frequencies  that  causes  the 

time  delayed  phase  shift  sensed  at  the  receiver  and  known  as 
the  Lloyd's  mirror  interference  phenomenon.  Since  both  the 
direct  and  reflected  rays  are  traveling  at  the  same  sound 
velocity  at  the  receiver,  the  ray  path  length  difference  can 
be  directly  converted  to  the  frequency  difference,  A  f  as 
observed  at  the  receiver  and  seen  in  the  Lloyd's  mirror 
interference  pattern. 

The  intensity  ratio  at  the  receiver  is  given  by 


/  ^ 

Intensity  ratio  ■  1+m2  |  — 1 


-  2  m 


COS  (  U)  r  ) 


t  N»7Vf 


The  ratio  reaches  a  minimum  value  corresponding  to 

maximum  interference  when  cos(wr)  =  1.  This  occurs  when 

«  t  »  2  »  or  any  integer  multiple  of  2r  . 

Note  that  X  =  c/f  where  X  =  wavelength  (m) 

c  a  phase  speed  (m/s) 
f  a  frequency  in  (hz) 

k  =  wave  number  (m-1 ) 


and  k  =  2  *•  / X 

Now  let  A0  equal  the  phase  difference  seen  at  the 
receiver  due  to  the  two  path  length  differences  and  defined 
as. 

At)  a  k  (SR-  SD),  where 
SR  a  Path  length  of  the  reflected  ray  and 
SD  >  Path  length  of  the  direct  ray. 

Let  AS  «  (SR-  SD) 

So  that  the  first  null  occurs  when 


A  0  a  k  (AS)  *  2  T  a  W  T 

Sok  (AS)  b«tb  2  t  =  2  *  (  A  S )  /  X  or 

A  S  «  2w\/2  w  3  X 

therefore  AS  =  c/f-j 

so  the  first  null  occurs  at  frequency 
f1  a  c/  AS 


The  second  null  occurs  at 
w  r  a  4  w  so  that 
AS  ■  4»X/2*  a  2  X  ■  2c/f 2 


therefore  second  null  occurs  at  frequency 


which  shows  that  f2  -  f-j  ■  f-j 
and  so  on.  So  that  for  the  n^1  null 
fn  »  nc/A  S  *  nf1 

and  fn  -  fn_i  «  f-|  or  more  generally 


where  m  >>  n  for  best  accuracy  so  that  the  measured 
difference  in  frequency  between  any  two  successive  nulls  can 
be  related  to  the  difference  in  path  length  required  to 
produce  that  null  by  the  equation 

_  C ( receiver ) 

f  a  — — 

AS 

where  C^receiverj  is  the  actual  sound  velocity  at 
the  receiver  at  the  time  of  reception. 


IV.  COMPARISON  of  the  ANALYTICAL  and  COMPUTER  SOLUTIONS  for 
the  ISOVELOCITY  SOUND  PROPAGATION  CASE 


This  analysis  demonstrates  the  accuracy  of  the  computer 
solution  by  using  an  isovelocity  sound  propagation  model. 
The  constant  sound  velocity  condition  results  in  straight 
line  transmission  paths  as  the  sound  travels  from  the  source 
to  the  receiver.  By  solving  the  geometric  model  analytically 
and  comparing  the  results  to  the  computer  output  one  can  see 
the  accuracy  of  the  Romberg  integration  technique  coupled 
with  the  four  integrals  derived  in  section  III  .  The 
isovelocity  case  is  the  simplest  propagation  model  yet  it 
exercises  all  phases  of  the  computer  program's  main 
calculation  algorithm.  The  only  section  of  the  overall 
program  not  tested  by  the  isovelocity  case  are  the  sound 
velocity  calculations  performed  as  the  user  inputs  a  varying 
tesiperature  profile.  Table  IV-1  compares  the  difference  in 
path  lengths  between  the  direct  and  reflected  sound  paths  as 
calculated  analytically  and  by  the  computer.  Of  the  20  data 
points  compared,  the  computer  missed  the  exact  ray  path 
length  for  a  horizontal  range  of  750  meters  by  only  19.16 
centimeters  in  the  worst  case  and  averaged  5.02  centimeters 
error  over  all  20  cases.  Table  IV-2  shows  the  resulting 
horizontal  range  error  due  to  the  error  of  the  computer's 


calculation  presented  in  Table  IV-1.  The  average  range  error 
for  the  20  points  is  only  1.52  meters. 

The  computer  program  steps  through  sea  depths  of  0  to 
200  meters  in  10  meter  increments  where  each  depth 
represents  a  possible  target  depth.  To  calculate  a  complete 
set  of  Snell  constants,  initial  sound  ray  angles  of  3  thru 
177  degrees  in  increments  of  3  degrees  (  measured  from  the 
vertical  )  are  used.  Ultimately  the  integration  section 
derives  a  path  length  difference  between  the  surface 
reflected  ray  and  the  direct  path  ray  for  each  of  these 
selected  starting  ray  direct  ions. The  59  initial  rays  times 
the  20  depth  increments  yields  a  total  of  1,180  data  points. 
These  point  solutions  are  used  to  interpolate  values  of  path 
length  difference  for  points  every  10  meters  of  horizontal 
range  out  to  1000  meters  which  results  in  a  two  dimensional 
matrix  containing  2000  (  20  depth  increments  times  100  range 
increments  )  horizontal  range  vs  depth  elements. 

To  compare  the  computer  solution  to  actual  geometric 
results,  five  ranges  for  each  of  four  depths  were  chosen  at 
random.  In  all  comparisons  the  receiver  is  set  at  a  depth  of 
27  meters  (approximately  90  feet). 
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TABLE  IV- 1 


ISOVELOCITY  CASE 
COMPARISON  BETWEEN 
GEOMETRIC  AND  COMPUTER  SOLUTIONS 


PATH  LENGTH  DIFFERENCE 
ANALYTICAL 


COMPUTER 


SOURCE 

DEPTHS  HORIZONTAL  RANGES.  METERS 

METERS 


10 


SO 


90 


130 


170 


_ LM _ _ Z5£ 


15.8107 

15.8228 

49.0755 

49.0196 

52.5987 

52.6316 

53.3450 

53.3808 

53.6208 

53.5714 


5.1908 

5.1939 

23.5992 

23.6220 

35.7218 

35.7143 

42.5840 

42.6136 

46.4312 

46.4396 


1.4170 

1.4418 

7.0274 

7.0389 

12.4171 

12.4688 

17.4439 

17.5644 

22.0132 

22.0264 


0.7195 

.6693 

3.5897 

3.5714 

6.4289 

6.2267 

9.2169 

9.0253 

11 .9302 
11.7739 
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TABLE  IV- 2 

COMPARISON  Of  PROGRAM  RANGES  vs.  ACTUAL  RANGES 


ACTUAL 

SOURCE 


HORIZONTAL 

RANGES. 

METERS 

20 

100 

380 

750 

10 

19.954 

99.935 

373.42 

806.3 

50 

20.134 

99.875 

379.37 

753.87 

90 

19.745 

100.038 

378.33 

774.84 

130 

19.42 

99.815 

377.075 

766.41 

170 

21  .26 

99.93 

379.725 

760.475 

COMPUTER 


HORIZONTAL  RANGES,  METERS 


COMPUTER  PROGRAM  FLOW  CHART 


A.  FLOW  CHART  and  PROGRAM  BACKGROUND 

The  following  block  diagram  flow  chart  demonstrates  the 
specific  steps  required  to  obtain  an  accurate  horizontal 
range  vs.  depth  solution  using  the  Lloyd  mirror  interference 
phenomenon.  Specific  values  or  data  used  in  the  included 
program  will  be  given  at  the  bottom  of  each  appropriate 
block. 

The  actual  numerical  integration  section  of  the  program 
is  taken  directly  from  the  Romberg  integration  program 
provided  by  Miller  [Ref.  7:  pp.  215-218]. 

The  sound  velocity  equation  used  in  the  program  was 
developed  by  Coppens  [Refs.  6  and  8]  and  selected  after 
several  equations  of  varying  complexity  were  considered. 
Specifically,  equations  derived  by  Wilson  and  Leroy  [Ref.  1: 
p.  105],  plus  Mackenzie  [Ref.  9]  were  compared  with  Coppens. 
All  results  were  within  1  to  2  centimeters  per  second  of 
each  other.  Coppen's  equation  was  chosen  because  of  the 
ability  to  correct  for  not  only  temperature,  pressure  and 
salinity,  but  also,  latitude  variations. 
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FLOW  CHAR? 


1 .  Initialization  and  Main  Menu 


.  INITIALIZATION 

1. 

DEFINE  ALL  VARIABLES 

2. 

DIMENSION  ALL  ARRAYS 

3. 

DEFINE  TOLERANCE 

DESIRED  FOR  NUMERICAL 

INTEGRATION 

I 


RANDOM  ACCESS  FILE 
CREATE  A  RANDOM 
ACCESS  FILE  TO 
STORE  THE  RANGE 
VS.  DEPTH  MATRIX 
SOLUTION 

I 


DISPLAY  MENU  OF 

PROGRAM  SUBROUTINES 
ALLOWING  USER  TO 

CHOOSE  DESIRED 
OPERATION  AS 

FOLLOWS: 

t. 

REINITIALIZE  RANDOM 
ACCESS  ANSWER  FILE 

AFTER  INITIAL 

RUN  AND  PRIOR  TO 

ALL  SUBSEQUENT  RUNS 

«•  • 

INPUT  WATER  TEMPERATURE 
INFO  AS  A  FUNCTION  OF 
DEPTH 

3. 

INPUT  SEA  SURFACE  TEMP¬ 
ERATURE  AND  AVERAGE 
GRADIENT 

4. 

INPUT  SOUND  VELOCITY 
PROFILE  AS  A  FUNCTION 

OF  DEPTH 

5. 

COMPUTE  HORIZONTAL  RANGE 
VS.  DEPTH  SOLUTION  AND 
STORE  IN  ANSWER  FILE 

6. 

DISPLAY  ANSWER  FILE  ON 

THE  CRT  SCREE:;  OR  LIST  ; 
OK  A  PRINTER 


7.  EXIT  PROGRAM  AND  RETURN 
TO  OPERATING  SYSTEM 


2.  Reinitialize  Random  Access  Answer  File 


9 


REINITIALIZATION 
DESIRED,  YES  OR  NO 


"wn  inni 

THIS  SECTION  WILL 

RETURN 

DESTROY  CURRENT 

NO* 

TO  MENU 

ANSWER  FILE 

QUERY  USER 

i 

k 

® 


-  -•  •  Jl 


Input  Water  Temperature. Profile 


INPUT  WATER  TEMPERATURE  IN 
DEGREES  CELSIUS  AT  UNIFORM 
INCREMENTS  OF  DEPTH  IN  METERS 
BEGINNING  WITH  SURFACE  (PROGRAM 
STEPS  IN  5  METER  INCREMENTS 
DOWN  TO  200  METERS) 
CALCULATE  SOUND  VELOCITY  IN  METERS 
PER  SECOND  FOR  EACH  1  METER 
INCREMENT  OF  DEPTH  AND  STORE 
ANSWER  IN  TEMPORARY  ARRAY 


Input  Sea  Surface  Temperature  And  Average  Gradient 


Input  Sound  Velocity  Profile 


6.  Compute  Range  vs.  Depth  Solution  And  Store  In  Random 
Access  File  ” 


(D 


INPUT  HYDROPi 
DEPTH  II 
FEET  OI 

IONE/ RECEIVER 

J  EITHER 
*  METERS 

>  —  - , 

LOOP  THROUGH  EACH  POSSIBLE 


[SOURCE  DEPTH,  (THIS  PROGRAM 
iUSES  10  METER  INCREMENTS  OF 
DEPTH  DOWN  TO  200  METERS) 

COMPUTING  PATH  LENGTH 
DIFFERENCES  FOR  EACH 
INCREMENT  OF  HORIZONTAL 
RANGE,  (THIS  PROGRAM  USES 
10  METER  INCREMENTS  OF 
HORIZONTAL  RANGE  OUT 
TO  1000  METERS) 

BY  DOING  THE  FOLLOWING: 

1  .  DETERMINE  THE  UPPER 

AND  LOWER  LIMITS  FOR  THE 
FOUR  SEPARATE  NUMERICAL 
INTEGRATIONS 

2.  LOOP  THROUGH  EACH 

POSSIBLE  INITIAL  RAY 

ANGLE  FROM  EACH  SOURCE 

DEPTH  COMPUTING  THE  FOLLOWING: 

a.  RETRIEVE  THE 

REQUIRED  SOUND  VELOCITIES 
FOR  EACH  INTEGRATION 

b.  COMPUTE  EACH 

UNIQUE  SNELL  CONSTANT 
GIVEN  AN  INITIAL  RAY 
ANGLE  AND  SOUND  VELOCITY 
FOR  THE  SOURCE  DEPTH 

C.  NUMERICALLY  SOLVE 

THE  FOUR  INTEGRALS  DERIVED 
IN  THE  MATHEMATICAL 
SOLUTION  DISCUSSION, 

L  —  —  tcont.  p.  9) - 
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cont.  from  p 


USING  A  ROMBERG  IMPROVE¬ 
MENT  [REF.  7]  ON  THE 
TRAPAZOI DAL  NUMERICAL 
INTEGRATION  TECHNIQUE 

d.  IF  INTEGRATION  ERROR 
OCCURS  WARN  USER  AND 
RETURN  TO  MENU 

e.  STORE  EACH  INTEGRATION 
ANSWER  IN  THE  APPROPRIATE 
TEMPORARY  ARRAY 

f.  CONTINUE  TO  LOOP  THROUGH 
ALL  POSSIBLE  RAY  ANGLES 
FROM  THE  SOURCE.  (THIS 
PROGRAM  USES  ANGLES  FROM 
3  DEGREES  FROM  THE  VERT¬ 
ICAL  TO  177  DEGREES  FROM 
THE  VERTICAL,  IN  3  DEGREE 
INCREMENTS) 

USE  LINEAR  INTERPOLATION  TO 
CALCULATE  THE  f ( DIFFERENCE) 
FOR  EACH  INCREMENT  OF  RANGE 
FOR  A  GIVEN  DEPTH 

STORE  ANSWERS  IN  RANDOM 
ACCESS  FILE 

CONTINUE  TO  LOOP  THROUGH 
DEPTH  INCREMENTS  UNTIL 
DONE 


RETURN 
TO  MENU 


I 


DISPLAY  OUTPUT  MENU  TO 
USER  ALLOWING  FOR  SELECTION 
OF: 

1.  CRT  SCREEN 

2 .  PRINTER 

RETRIEVE  HORIZONTAL  RANGE  VS 
DEPTH  ANSWER  MATRIX  FROM 
RANDOM  ACCESS  FILE  AND 
OUTPUT  TO  USERS  CHOICE 
FROM  ABOVE 


RETURN 
TO  MENU 


<s> 
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WARN  USER  THAT  THIS 
OPERATION  WILL  DESTROY 
ALL  TEMPORARY  FILES 
SUCH  AS  TEMPERATURE 
AND  SOUND  VELOCITY 
PROFILES 


ASK  USER  IF 

EXIT  IS  STILL 

NO  t 

RETURN 

DESIRED, 

YES  OR  NO 

TO  MENU 

YES 

1 

4 

EXIT  TO  OPERATING 


SYSTEM 


TT-Tr^TTr*  .*  »~:r-  vr tttt 

1 

VI. 

DEFINITION  Of  COMPUTER  PROGRAM  VARIABLES 

s 

l 

ADO 

Equation  used  to  compute  raypath  length  in 

meters  for  both  the  direct  and  reflected  paths. 

AFLAC 

Branching  flag  used  in  solving  for  answer  matrix 

elements 

AMS 

Answer  obtained  from  numeric  integration 

4 

AOCJT 

Answer  matrix  element  output  to  CRT  screen  or 

printer 

BFLAG 

Branching  flag  used  in  solving  for  answer  matrix 

t 

elements 

t 

BNDLOW 

Lower  limit  for  evaluating  numeric  integration 

* 

II 

\ 

\ 

BOUNDHI 

Upper  limit  for  evaluating  numeric  integration 

| 

BODNDLO 

Direct  path  lower  limit 

' 

CLAT 

Latitude  variation  used  in  sound  velocity 

. 

• 

equation 

> 

> 

_5 

DDEEP 

First  half  of  depth  effect  used  in  sound 
velocity  equation 

DEEP 

Second  half  of  depth  effect  used  in  sound 
velocity  equation 

DEG. 

Average  whole  degree  latitude  of  operating  area 
used  to  calculate  CLAT 

DELTA 

Line  integral  uniform  divisions  used  in  Romberg 

c 

numeric  integration  routine 

£ 

DHZ 

Direct  path  horizontal  range  upper  bound  used  in 

* 

interpolating  direct  path  length  for  desired 

> 

range  increment 

.j 

DIRFRAC 

Direct  path  interpolation  fraction 

-I 

l 

3 

DIRRAY 

Direct  path  interpolated  path  length 

a 

\ 

> 

i' 

• 
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DLONG 


DLOW 

DPH 

DSHORT 

F ANSWER 

FINAL  ( 

GRAD 

HYDEEP 

HYDPTH 

HYSPD 

I 

IFLAG 

II 
IN 

INCR 

INN 

J 

JMP 

K 

N 

METER 


Dirct  path  length  upper  bound  used  in 
Interpolating  path  length  for  desired  range 
increment 

Direct  path  hoirzontal  range  lower  bound  used  in 
interpolating  direct  path  length  for  desired 
range  increment 

Depth  increment  used  for  outputing  answer  matrix 

Direct  path  length  lower  bound  used  in 
interpolating  path  length  for  desired  range 
increment 

FL1  Random  access  file  used  to  store  the  horizontal 
range  vs.  depth  answer  matrix 

)  Answer  array  used  to  temporarily  hold  the 
interpolated  values  before  storing  them  in 
random  access  file 

Average  temperature  gradilent  input  by  user 
Hydrophone/receiver  depth  in  feet 
Hydrophone /receiver  depth  in  meters 
Sound  velocity  at  hydrophone/receiver  depth 
Loop  increment 

Loop  increment  used  to  branch  through  one  of 
four  integrated  equations 

Loop  increment 

Loop  increment 

Loop  increment 

Loop  increment 

Counter 

Array  element  pointer 

Counter 

Counter 

Array  element  pointer 


vT  *7*  <m  r  v* 


v*. 


j 


MSTEP 

Loop  increment 

MSV 

Array  element  pointer 

MTRHI 

Array  element  pointer 

MTRLO 

Array  element  pointer 

MTRS 

Array  element  pointer 

MTRSTEP 

Array  element  pointer 

MULT 

Random  access  file  pointer 

N 

Array  element  pointer 

NANGL 

Ray  angle  from  source 

NMTR 

Horizontal  range  increment 

MN 

Counter 

NSTEP 

Ray  angle  from  source 

NT 

Array  element  pointer 

NTRA 

Array  element  pointer 

SUM 

Main  menu  index  selection 

NUMB 

Output  menu  index  selection 

NUPLMT 

Incremented  source  depth 

NX  (  ) 

Array  of  counters  used  during  integration 
accuracy  testing  loop 

PART 

Partitioning  factor  used  in  numerical 
integration  routine 

PIECES 

Incremented  divisor  used  to  partition  depth 
integration  limits 

RAD 

Conversion  of  ray  angle  from  degrees  to  radians 

RAYDIR  ( 

)  Direct  ray  path  length  calculated  from  equation 
#12  page  25 

RAYREF 

Reflected  path  length  calculated  from  equation 
#13  page  25 
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REFFRAC  Reflected  path  interpolation  fraction 

REFLO  Reflected  path  lower  limit.  Negative  value  of 

bound lo 

REFRAY  Reflected  path  interpolated  path  length 

RHI  Reflected  path  horizontal  range  upper  bound  used 

in  interpolating  reflected  path  length  for 
desired  range 

RLONG  Reflected  path  length  upper  bound  used  in 

interpolating  path  length  for  desired  range 
increment 

RLOW  Reflected  path  horizontal  range  lower  bound  used 

in  interpolating  reflected  path  length  for 
desired  range 

RMTR  Depth  increment  for  interpolation  routine 

RNGDIR  (  )  Direct  path  horizontal  range  array  calculated 
from  equation  #7  page  21 

RNGREF  (  )  Reflected  path  horizontal  range  array  calculated 
from  equation  §9  page  22 

RSHQRT  Reflected  path  length  lower  bound  used  in 

interpolating  path  length 

SALIN  Salinity  value  for  operating  area  input  by  user 

SNDHI  Sound  velocity  at  numeric  integral  upper  bound 

SNDLO  Sound  velocity  at  numeric  integral  lower  bound 

SNDSPS  (  )  Sound  velocity  profile  array  input  directly  by 
user  or  calculated  from  temperature,  salinity, 
and  latitude  inputs 

SNELL  Snell  constant  equation  result  calculated  from 

initial  ray  angle  and  depth  of  source 

SOUND  Sound  velocity  at  incremented  source  depth 

SST  Sea  surface  temperature 

START  Initial  estimate  of  integral  answer 

Equation  used  to  compute  horizontal  range  in 
meters  for  both  the  direct  and  reflected  paths 


SUM 


SVP 
TEMP 
TEST  (  ) 

TGTDPTH 

TOL 

TOT 

TTEMP 

WALK 

X 


Input  sound  velocity  profile  by  user 

Degrees  Celsius  temperature  divided  by  10 

Successive  numeric  integration  answer  array  used 
to  check  accuracy  of  answer 

Incremented  source  depth 

Accuracy  desired  for  numerical  integration 
routine  answer 

Summation  of  partitioned  numerical  ilntegration 
values 

Input  temperature  used  in  calculating  sound 
velocity  profile 

Integration  routine  smoothing  factor 
Incremented  depth  value  for  each  partition 
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APPENDIX  A 


Urlck  [Ref.  1:  p.  124]  briefly  discusses  the  Intensity  of  sound 
reaching  a  receiver  via  a  direct  path  and  a  surface  reflected  path.  The 
following  derives  the  ratio  of  the  total  Intensity  of  the  signal  received 
at  a  point  when  multipath  Lloyd's  mirror  interference  is  present,  to  the 
total  Intensity  resulting  only  from  the  direct  path. 


<I>  2  /  R0\2  /*d\ 

‘ 1  *  “  (kj)  -  cos(“T ) 


WHERE:  <I>,  the  time  averaged  Intensity  and  I,  the  Instantaneous 
Intensity  are  defined  as  follows: 


(1) 

(2) 


0 


(P0(r)  ♦  PR(r))2 
- - Zpc - 


T  ■  1  period 


PQ(r)  and  PR(r)  are  both  functions  of  time 
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GIVEN: 


W  Deceived  *  pdirect  +  Reflected 

p0  up0 

*  k—  COS  Hit  -  75 -  COS 

kd  kr 

where:  u  ■  surface  pressure  reflection  coefficient 

t  *  time  delay  caused  by  the  difference  in  the  two  path 
lengths 

2  P^ 

(2)  ‘•■fe-i* 

where:  P  »  RMS  pressure  amplitude  of  the  direct  path  ray 
NOTE:  It  Is  assumed  that  the  transmission  loss  is  due  only  to  spherical 


spreading  and  that  absorption  in  the  medium  is  negligible. 


SOLUTION: 


Pg  -wPq 

PDIR  *  C0S  {wt)  PREF  *  COsCcu(t-r)3 


PRECEIVE0  ’  P0IRECT  +  PREFLECTED 


RECEIVED 


*  Po[^“  COS(  Hit)  -  COS  ( (jj(t-x)  )j 


p£  ,  ^0  cos^fiiit 
pC  2pC  ^2 


♦  cos  ( w(t-T)) 

kr 


-  ‘  COS  (u>t)  COS  Do(t-x)] 

krko 
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•lm.  :  »'v- 


SINCE:  <I> 


♦/■ « 


AND  u  *  2*f  *  T  *  — 

T  u 

2k  fm 

<!>-#■/  I  dt 


*/ 


<I> 


U>  P 


2wpc 


2w/u» 


cos^(ut)  dt 


2  ,2,/“ 
*%/  ' 


cos ^(w(t-T ) )  dt 


2k/ a) 

cos(wt)  C0s(u(t”T ) )  dt 

0 

Solving  the  three  Integrals  separately 


2*/b> 

I" 


cos  volt)  dt 


h  +  is s1n(2 

% 


2k/  “ 
0 


sin(4ir)] 


“«D 
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Assembling  the  answers  to  the  three  Integrals  yields: 


<I> 


“>R5 

2*pC 

r>2  r 


+  ^  cos<“rt 


o»RJ  “''R''D 


<i> 


2pcRf 


1  + 


M  2yRr 


RR  R 


COS(  (jx) 


NOTE: 


<I0>  =* 


2pCr; 


<I>  «  <I0> 


Finally 


2b2 


u  R5  2mR0 

1  +  — -  ~k —  cos(ut) 


This  result  shows  the  dependence  of  the  Intensity  at  the  receiver  on 
the  difference  In  path  lengths  RqIRECT  and  rreplecjed  and  the  time  delay, 
r,  caused  by  this  difference.  This  Intensity  ratio  Is  a  minimum  (corre¬ 
sponding  to  a  null  on  the  Lloyd's  mirror  interference  pattern)  when 
cos  (i»t  *  1  or  (at  *  Integer  multiples  of  2w. 


COMPUTER  OUTPUT 


A  sample  of  the  actual  computer  output  is  listed  below. 
Data  values  are  frequency  measurements  Af  (Hertz)  between 
successive  nulls  on  the  Lloyd*  s  mirror  pattern  (see  Figure 
II-l). 

Given  a  particular  source  depth  and  a  Af,  the  source  range 
can  be  determined  from  the  corresponding  output  column.  The 
source  depth  can  be  determined  from  the  appropriate  output  row 
if  the  horizontal  range  to  the  source  from  the  receiver  and 
the  Af  are  known. 

Group  of  Ranges  From  10  to  70  Meters 
Depth  of  10  Meters 


80.6 

94.8  114.1 

136.0 

159.8 

184.4 

210.1 

Depth  of 
41.9 

20  Meters 

50.2  60.1 

71.0 

82.7 

94.8 

107.3 

Depth  of 
31.6 

30  Meters 

37.3  43.8 

50.8 

58.3 

66.1 

74.3 

Depth  of 
29.2 

40  Meters 

32.5  36.8 

41.7 

46.9 

52.5 

58.3 

Depth  of 
28.5 

50  Meters 

30.6  33.5 

36.9 

40.8 

44.9 

49.4 

COMPUTER  PROGRAM 


10  DEFINT  I,J,K,L,M,N 

20  DIM  NX(16),TEST(136),SNDSPD(201),RNGDIR(180) 

30  DIM  RNGREF(180),RAYDIR(180),RAYREF(180),FINAL(1  01  ) 

40  TOL  »  .05 

50  OPEN  "R",#1  ,"FANSWER.FL1",4 
60  FIELD  #1  ,  4  AS  FI  $ 

100  PRINT  CHR$( 26):  PRINT  "SELECT  OPERATION": PRINT 

110  PRINT  1  , "INITIALIZE  ANSWER  FILE, DO  AT  THE  BEGINNING  OF 

EACH  NEW  BT  INPUT" 

120  PRINT  2, "INPUT  BT  INFO  IN  5  METER  INCREMENTS" 

130  PRINT  3, "INPUT  SEA  SURFACE  TEMP  AND  GRADIENT  FOR  FIRST 
200  METERS" 

135  PRINT  4f "INPUT  SOUND  VELOCITY  PROFILE  IN  5  METER 
INCREMENTS" 

140  PRINT  5, "CALCULATE  LLOYD'S  MIRROR  HORIZONTAL  RANGE  .VS. 
DEPTH  CHART" 

150  print  6, "Output  chart  to  printer" 

160  PRINT  7,  "EXIT  PROGRAM" 

170  PRINT:PRINT:INPUT  "ENTER  THE  NUMBER  OF  THE  OPERATION 
DESIRED  ",NUM 

180  IF  (NUM  <  1  )  OR  (NUM  >7 )  THEN  PRINT  "IMPROPER  INPUT":GOTO 
100 

190  ON  NUM  GOSUB  1000,2000,5000,5500,6000,9600,10000 
200  GOTO  100 

1000  REM  INITIALIZE  ANSWER  FILE 

1010  PRINT  CHR${ 26) SPRINT  "ARE  YOU  SURE  YOU  WANT  TO  ZERO  OUT 
THE  ANSWER  FILE" 

1030  INPUT  "ENTER  Y  FOR  YES,  N  FOR  NO  ",A$ 

1040  IF  A$<  >"Y"  THEN  RETURN 

1050  CLOSE  #1 

1060  KILL  "FANSWER.FL1" 

1070  OPEN  "R",#1  , "FANSWER.FL1", 4 
1080  FIELD  #1,4  AS  FI  $ 

1090  RETURN 

2000  REM  INPUT  BT  INFO  IN  5  METER  INCREMENTS, THEN  CALCULATE 
SOUND  VELOCITY 

2010  REM  USING  ALAN  R.  COPPENS  EQUATION  PRESENTED  IN  J. 
ACOUS.  SOC.  OF  AMER. 

2015  REM  VOLUME  69,  PAGE  862,  1981  AND  PUBLISHED  IN 

FUNDAMENTALS  OF  ACOUSTICS 

2016  REM  THIRD  EDITION,  BY  KINSLER  FREY  COPPENS  AND 
SANDERS 

2020  GOSUB  3010:  GOSUB  4000:PRINT  CHR$(26) 

2030  PRINT  "ENTER  BT  INFO  IN  5  METER  INCREMENTS  STARTING 
WITH  THE  SURFACE" 

2040  FOR  I  «  0  TO  200  STEP  5 

2050  PRINT:PRINT  I  "METERS": INPUT  "TEMP  IN  CELSIUS  » 
",TTEMP: PRINT 

2060  PRINT  "YOU  INPUT  "TTEMP"  DEGREES  CELSIUS": 


INPUT  "IS  THIS  CORRECT  7  ENTER  Y  FOR  YES  N  FOR  NO 

";b$ 

2070  IF  B$<  >"Y"  THEN  PRINT  "REENTER  TEMP":GOTO  2050 
2075  TEMP  ■  TTEMP  /  10! 

2080  METER  *  I 

2100  DEEP  »  (16.23  +.253*TEMP)*  CLAT  *  CSNG(I)  /1 000!  +(.213- 

.1*TEMP)  * 

(CLAT  *  (CSNG(I)  /1000l))*2 

2105  DDEEP  ■  (.01  6  +  .0002  *  (SALIN-35))  *  ( S ALIN - 

35 ) *TEMP*CLAT*CSNG ( I ) /I  000 l 

2110  SNDSPD(METER) -1449.05  +  45.7  *  TEMP  -  5.21  *  (TEMPa2 )  + 
^2  3*  (TEMP* 3)  + 

(1 .333-.126*TEMP+.009*(TEMP*2) )*(SALIN- 
35 ) + ( DEEP+DDEEP ) 

2111  IF  I  «  0  GOTO  2120 

2112  SLOW  »  SN DS  PD  (  I  -  5  ) 

2113  SHI  -  SNDSPD(  I  ) 

2114  FOR  J  »  1  TO  4 

2115  SNDSPD(  I  -  5  +  J  )  =  SLOW  +  (  J  /  5  )  *  (  SHI  -  SLOW  ) 

2116  NEXT  J 
2120  NEXT  I 
2130  RETURN 

3000  REM  ENTRY  OF  OPAREA  LATITUDE 

3010  PRINT  CHR$(  26):  PRINT  "ENTER  AVERAGE  DEGREES  LATITUDE  OF 
OPAREA" 

3020  INPUT  "IN  WHOLE  DEGREES  ",DEG 

3030  PRINTsPRINT  "YOU  ENTERED  "DEG"  DEGREES ,  IS  THIS  CORRECT 
7" 

3040  INPUT  "ENTER  Y  FOR  YES  N  FOR  NO  ",CC$ 

3050  IF  CC$<>"Y"  THEN  PRINT  "REENTER  LATITUDE  ":GOTO  3010 
3060  CLAT  -  1  -  .0026  *  (COS(  DEG  *  3.14159  /  901  )) 

3070  RETURN 

4000  REM  ENTRY  OF  SALINITY  VALUE  FOR  SOUND  SPEED 

COMPUTATIONS 

4010  PRINT  CHR$(26) 

4020  PRINT  "ENTER  THE  SALINITY  VALUE  AS  PARTS  PER  THOUSAND 
FOR  THE  OPAREA" 

4030  INPUT  "IF  UNKNOWN,  ENTER  A  VALUE  OF  35,  SALINITY  * 
" , SALIN 

4040  PRINTsPRINT  "THE  SALINITY  LEVEL  IS  "SALIN"  PARTS  PER 
THOUSAND" 

4050  INPUT  "IS  THIS  CORRECT  7  ENTER  Y  FOR  YES  N  FOR  NO 

",bb$ 

4060  IF  BB$  <>  "Y"  GOTO  4010 
4070  RETURN 

5000  REM  ENTER  SEA  SURFACE  TEMP  AND  AVERAGE  GRADIENT  FOR 
FIRST  200  METERS 
5010  PRINT  CHR$(26) 

5020  INPUT  "ENTER  THE  SEA  SURFACE  TEMPERATURE  IN  DEGREES 
CELSIUS  ",SST 

5030  PRINTSPRINT  "YOU  ENTERED  "SST"  DEGREES  CELSIUS  FOR  THE 


,  .v.'  ,\v 


••  J, 
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SEA  SURFACE  TEMP" 

5040  INPUT  "IS  THIS  CORRECT  ?  ENTER  Y  FOR  YES  N  FOR  NO 
\C$ 

5050  IF  C$  <>  "Y"  THEN  PRINT  "REENTER  SEA  SURFACE  TEMP": GOTO 
5020 

5060  PRINT 

5070  PRINT  "ENTER  THE  TEMPERATURE  GRADIENT  (DEG 
CELSIUS/METER) " 

5075  INPUT  "FOR  THE  FIRST  200  METERS  ",  GRAD:  PRINT 

5080  PRINT  "YOU  ENTERED  A  TEMPERATURE  GRADIENT  OF  "GRAD" 

DEGREES  CELSIUS/METER" 

50  90  INPUT  "IS  THIS  CORRECT  ?  ENTER  Y  FOR  YES  N  FOR  NO 
"  ,D$ 

5100  IF  D$  <>  "Y"  THEN  PRINT  "REENTER  GRADIENT" :GOTO  5070 
5110  GOSUB  3010 
5120  GOSUB  4000 

5125  PRINT  CHR$(  26 ): PRINT  "PLEASE  WAIT,  CALCULATING  SOUND 
VELOCITY  VALUES" 

5130  FOR  I  »  0  TO  200 

5140  TEMP  «  (SST  +  (CSNG(I)  *  GRAD))  /  10! 

5150  DEEP  «  (16.23  +  .253*TEMP)*  CLAT  *  CSNG(I)  /1000!  + 

( .21 3-. 1  *TEMP) 

*  (CLAT  *  (CSNG(I)  /1 000!  ))a2 

5155  DDEEP  «  (.016  +  .0002  *  (SALIN  -  35))  *  (SALIN-35)  * 
TEMP  * CLAT  * 

CSNG(I)  /  1000! 

5170  SNDSPD(I)  »  1449.05  ♦  45.7  *  TEMP  -  5.21  *  (TEMPa  2 )  ♦ 
.23  *  (TEMP* 3  ) 

♦  (  1.333  -  .126  *  TEMP  +  ,009*TEMPa2  )*(SALIN  - 

35 ) +DEEP+DDEEP 
5180  NEXT  I 
5190  RETURN 

551 0  REM  ENTRY  OF  SOUND  VELOCITY  PROFILE 
5525  PRINT  CHR$(26) 

5530  PRINT  "ENTER  SOUND  VELOCITY  PROFILE  IN  5  METER 
INCREMENTS  " 

5540  PRINT  "STARTING  WITH  THE  SURFACE  SOUND  VELOCITY  IN 
METERS/SECOND" 

5550  FOR  II  *  0  TO  200  STEP  5 

5560  PRINTsPRINT  II  "METERS"  :INPUT  "  SOUND 
VELOCITY  (METERS/ SEC) -  ",SVP 

5570  PRINT: PRINT  "YOU  INPUT  "SVP"  METERS/SEC" 

5580  INPUT  "IS  THIS  CORRECT  ?  ENTER  Y  FOR  YES  N  FOR  NO 
";SS$ 

5590  IF  SS$<  >"Y"  THEN  PRINT  "REENTER  SOUND  VELOCITY":GOTO 
5560 

5610  SNDSPD(II)  *  SVP 

5611  IF  II  »  0  GOTO  5620 

5612  SLOW  -  SNDSPD(  II  -  5  ) 

5613  SHI  *  SNDSPD(  II  ) 

5614  FOR  JJ  -  1  TO  4 
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5615  SNDSPD(  II  -  5  ♦  JJ  )  -  SLOW  ♦  (  JJ  /  5  )  *  (  SHI  - 
SLOW  ) 

5616  NEXT  JJ 
5620  NEXT  II 
5630  RETURN 

6000  REM  CALCULATE  LLOYD'S  MIRROR  HORIZONTAL  RANGE  .VS. 
DEPTH  CHART 

6005  PRINTsPRINTsPRINT  "WILL  HYDROPHONE  DEPTH  BE  INPUT  IN 
FEET  OR  METERS  ?" 

6006  INPUT  "ENTER  F  FOR  FEET  OR  M  FOR  METERS  ",  AHYD$ 

6007  IF  AHYD$  »  "F"  GOTO  6010  ELSE  IF  AHYD$  »  "M"  GOTO  6  052 

6008  PRINT: PRINT: PRINT  "INPUT  ERROR"  :GOTO  6005 

6010  PRINT  CHR$(  26 ): INPUT  "ENTER  HYDROPHONE  DEPTH  IN  FEET 
w  HYDEEP 

6020  PRINT :PRINT  "YOU  ENTERED  A  HYDROPHONE  DEPTH  OF  "HYDEEP" 
feet" 

6030  INPUT  "IS  THIS  CORRECT  ?  ENTER  Y  FOR  YES  N  FOR  NO 

"  ,AA$ 

6040  IF  AA$  <>  "Y"  GOTO  6010 

6050  HYDPTH  -  (HYDEEP  *  .3048) 

6051  GOTO  6060 

6052  PRINT  CHR$(  26):  INPUT  "ENTER  HYDROPHONE  DEPTH  IN  METERS 
", HYDPTH 

6053  PRINT:PRINT  "YOU  ENTERED  A  HYDROPHONE  DEPTH  OF 
"HYDPTH"  METERS" 

6054  INPUT  "  IS  THIS  CORRECT  ?  ENTER  Y  FOR  YES  N  FOR  NO 

",BHY$ 

6055  IF  BHY$  <>  "Y"  GOTO  6052 
6060  HYSPD  -  SNDSPD(  HYDPTH  ) 

61  70  FOR  NUPLMT  -  10  TO  200  STEP  1  0 
6175  PRINT  "NUPLMT-"  NUPLMT 
6180  TGTDPTH  -  CSNG(  NUPLMT  ) 

6185  IF  HYDPTH  «  TGTDPTH  THEN  HYDPTH  «  HYDPTH  -  3 

6190  IF  TGTDPTH  >  HYDPTH  GOTO  6230 

6200  BOUNDLO  «  TGTDPTH 

6210  BOUHDHI  «  HYDPTH 

6220  GOTO  6250 

6230  BOUNDLO  »  HYDPTH 

6240  BOUNDHI  «  TGTDPTH 

6250  REFLO  «  -(BOUNDLO) 

6360  FOR  NSTEP  «  3  TO  177  STEP  3 

6365  PRINT  'ttSTEP-"  NSTEP 

6370  RAD  -  CSNG(NSTEP)  *  3.1  41  59  /  1  801 

6380  MTRSTEP  -  TGTDPTH 

6390  SOUND  *  SNDSPD(MTRSTEP) 

6400  MTRHI  -  BOUNDHI 
6410  SNDHI  -  SNDSPD(MTRHI) 

6420  MTRLO  -  BOUNDLO 
6430  SNDLO  -  SNDSPD(MTRLO) 

6440  SNELL  -  SIN(  RAD  )  /  SOUND 
6480  FOR  IFLAG  -  1  TO  4 
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7000  REM  ROMBERG  NUMERICAL  INTEGRATION  SUBROUTINE 
7015  START  ■  01 
7020  PIECES  -  11 
7030  NX(1)  -  1 

BNDLOW1  BOUIHJLO^"  *  ^  ^  0R  ^IFLAG*4^  then  BNDLOWsREFLO  else 

7041  IP  (I FLAG  »  2)  OR  (I FLAG  =  4)  GOTO  7045 

7042  IF  (TGTDPTH  >  HYDPTH)  AND  (NSTEP  >  90)  GOTO  8105 

7043  IF  (TGTDPTH  <  HYDPTH)  AND  (NSTEP  <  90)  GOTO  8105 

7044  GOTO  7050 

7045  IF  NSTEP  >  90  GOTO  8105 

7050  DELTA  -  (BOUNDHI  -  BNDLOW)  /  PIECES 

7051  DENOM1  •  (SNELL) *2  *  (SNDLO) *2 

7052  DENOM2  »  (SNELL)8 2  *  (SNDHI)a2 

DENOM2F»DE99941  =  1  1  raEN  DEN0M1  *  *999  :IP  DENOM2  =  1!  THEN 

7054  IF  (  DENOM1  >  11  )  OR  (  DENOM2  >  11  )  GOTO  8105 
7060  IF  (IFLAGsl)  OR  (IFLAG-2)  GOTO  7085 
7065  FOR  IN  «  1  TO  2 

7068  IF  IN  »  1  THEN  DENOM  «  DENOM1  ELSE  DENOM  *  DENOM2 
7070  START  -  START  ♦  ((11)/  SQR(1  1  -  DENOM)) 

7075  NEXT  IN 

7080  GOTO  7100 

7085  FOR  INN  «  1  TO  2 

7088  IF  INN  «  1  THEN  DENOM  »  DENOM1  ELSE  DENOM  »  DENOM 2 

7089  IF  INN  «  1  THEN  SPEED  «  SNDLO  ELSE  SPEED  *  SNDHI 

7090  START  «  START  ♦  ((SNELL  *  SPEED)  /  SQR(1  1  -  DENOM)) 

7095  NEXT  INN  ' ' 

71  00  TOT  *  START  /  2 1 

7110  TEST(T)  -  DELTA  *  START  /  21 
7120  N  «  1 
7130  NN  -  2 
7140  N  *  N  +  1 


7150  WALK  *  41 

7160  NX(N)  «  NN 

7170  MSTEP  »  CINT(  PIECES  ) 

7180  PIECES  «  PIECES  *  21 

7190  DELTA  «  (BOUNDHI  -  BNDLOW)  /  PIECES 

7500  REM  COMPUTE  TRAPEZOIDAL  SUM 

7510  FOR  I  NCR  ■  1  TO  MSTEP 

7520  PART  »  CSNG(INCR)  *  21  -  11 

7530  X  »  BNDLOW  ♦  (DELTA  *  PART) 

7540  MTRS  »  CINT  (  ABS  (X)  ) 

7545  IF  MTRS  >  200  THEN  MTRS  >  200 
7550  SPEED  «  SNDSPD(  MTRS  ) 

7552  DEMO  M3  ■  (SNELL) a2  *  (SPEED) a2 

7554  IF  DENOM 3  -  11  THEN  DENOM 3  «  .999 

7555  IF  DENOM 3  >  11  GOTO  8105 

7570  IF  (IFLAG>1)  OR  (IFLAG«2)  GOTO  7600 
7580  TOT  -  TOT  +  (  (11)  /  SQR(11  -  DENOM3)) 
7590  GOTO  7610 
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7600  TOT  »  TOT  +  (  (SNELL  *  SPEED)  /  SQR(1i  -  DENOM3 ) ) 

7610  NEXT  INCH 

7620  TEST(  NN  )  «  TOT  *  DELTA 
7630  K  •  N  -  1 
7640  NTRA  ■  NX(  K  ) 

7650  FOR  M  -  1  TO  K 
7660  J  «  NN  +  M 
7670  NT  «  NTRA  +  M  -  1 

7680  TEST(J)  *  (  TEST(J-I)  *  WALK  -  TEST(NT) )  /  (WALK-1 1) 
7690  WALK  -  WALK  *  41 
7700  NEXT  M 

7705  IF  (I FLAG  «  1  OR  IFLAG  =  2)  AND  (TEST(J)  >  30001)  GOTO 
8105 

7710  IF  N  <  5  GOTO  7760 

7720  IF  TEST(  NN  +  1  )  *  01  GOTO  7750 

7730  IF  ABS(TEST(NTRA+1  )-TEST(NN  ♦  1))  <*  (ABS(TEST(NN  + 
1)*TOL))  GOTO  8000 

7740  IF  ABS(  TEST(NN-I)  -  TEST(J))  <«  (ABS(  TEST(J)  *  TOL)) 
GOTO  8000 

7750  IF  N  >  15  GOTO  9500 
7760  NN  *  J  +  1 
7770  GOTO  7140 
8000  ANS  a  TEST(J) 

8030  ON  IFLAG  GOTO  8040,8060,8080,8100 

8040  RNGDIR (NSTEP )  a  ANS 

8050  GOTO  8105 

8060  RNGREF (NSTEP)  *  ANS 

8070  GOTO  8105 

8080  RAYDIR (NSTEP)  «  ANS 

8090  GOTO  8105 

8100  RAYREF (NSTEP)  a  ANS 

8105  NEXT  IFLAG 

8110  NEXT  NSTEP 

8115  REM  USE  LINEAR  INTERPOLATION  TO  CALCULATE  ANSWER  CHART 
VALUES  IN  RANGE 

8117  REM  AND  DEPTH  INCREMENTS  OF  10  METERS 

8120  FOR  NMTR  -  1  TO  100 

8130  RMTR  «  CSNG(  NMTR  )  *  101 

8135  DIRRAY  a  01  :  REFRAY  a  01 

8140  AFLAG  a  01 

8150  BFLAG  -  01 

8160  IF  TGTDPTH  <  HYDPTH  GOTO  8179 

8161  FOR  NANGL  a  3  TO  87  STEP  3 

8162  IF  RNGDIR(  NANGL  )  <«  RMTR  GOTO  8171 

8163  DLOW  a  RNGDIR  (  NANGL  -  3  ):DHI  -  RNGDIR(  NANGL  ) 

8164  IF  DLOW  a  DHI  THEN  DIRFRAC  -  0  .-GOTO  8166 

8165  DIRFRAC  -  (RMTR  -  DLOW)  /  (DHI  -  DLOW) 

81  66  IF  DLOW  >  0  GOTO  81  68 

8167  DSHORT  -  BOUNDHI  -  BOUNDLO  .-GOTO  8169 

8168  DSHORT  a  RAYDIR (  NANGL  -  3  ) 

8169  DLONG  a  RAYDIR(  NANGL  ) 
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81  70  DIRRAY  >  OSHORT  +  (  DLONG  -  D SHORT  )  *  DIRFRAC  :  GOTO 
8210 

8171  IF  (  NANGL  «  87  )  AND  (  RNGDIR(  87  )  >  0  )  GOTO  8172 
ELSE  GOTO  8177 

8172  DLOW  *  RNGDIR(  87  ) 

8173  IF  DLOW  -  RMTR  THEN  DIRRAY  *  RAYDIR(  87  )  :  GOTO  8210 

8174  DIRFRAC  «  (  RMTR  -  DLOW  )  /  (  1000  -  DLOW  ) 

8175  DSHQRT  -  RAYDIR(  NANGL  )  :  CLONG  =  SQR(  1000 a2  + 
(BOUNDHI  -  BOUNDLO ) a  2  ) 

81  76  DIRRAY  *  DSHORT  +  (  CLONG  -  DSHQRT  )  *  DIRFRAC  :  GOTO 
8210 

8177  NEXT  NANGL 

8178  GOTO  8410 

8179  FOR  MANGL  -  93  TO  177  STEP  3 

8180  IF  RMTR  >«  RNGDIR(  MANGL  )  GOTO  8181  ELSE  GOTO  8195 

8181  IF  (  MANGL  =  93  )  AND  (  RNGD IR(  93  )  >  0  )  GOTO  8182 
ELSE  GOTO  8186 

8182  DLOW  *  RNGDIR(  93  ) 

8183  IF  DLOW  -  RMTR  THEN  DIRRAY  *  RAYDIR(  93  )  :  GOTO  8210 

8184  DIRFRAC  *  (  RMTR  -  DLOW  )  /  (  1000  -  DLOW  )  :  DSHORT  * 
RAYDIR(  93  )  : 

BLONG  «  SQR  (  1  000a2  +  (  BOUNDHI  -  BOUNDLO  )a2  ) 

8185  DIRRAY*  DSHORT  +  (  BLONG  -  DSHORT  )  *  DIRFRAC  :  GOTO 
8210 

8186  DLOW  *  RNGDIR(  MANGL  )  :DHI  *  RNGDIR(  MANGL  -  3  ) 

8187  IF  (  DLOW  <*  RMTR  )  AND  (  DHI  >-  RMTR  )  GOTO  8188  ELSE 
GOTO  8195 

8188  IF  DLOW  »  DHI  THEN  DIRFRAC  *  0  :  GOTO  8190 
81  89  DIRFRAC  *  (  RMTR  -  DLOW  )  /  (  DHI  -  D^OW  ) 

81  90  IF  DLOW  >  0  GOTO  81  92 

8191  DSHORT  *  BOUNDHI  -  BOUNDLO  :  GOTO  8193 

8192  DSHORT  -  RAYDIR(  MANGL  ):  DLONG  *  RAYDIR(  MANGL  -  3) 

8194  DIRRAY  *  DSHORT  +  (  DLONG  -  DSHORT  )  *  DIRFRAC  :  GOTO 
8210 

8195  IF  (  MANGL  <>  177  )  GOTO  8199 

8196  DLOW  *  0  :  DHI  *  RNGDIR(  1  77  )  S  DIRFRAC  *  RMTR  /  DHI 

8197  DSHORT  *  BOUNDHI  -  BOUNDLO  :  DLONG  *  RAYDIR(  177  ) 

8198  DIRRAY  *  DSHORT  *  (  DLONG  -  DSHORT  )  *  DIRFRAC  :  GOTO 
8210 

8199  NEXT  MANGL 

8200  GOTO  8410 

8205  IF  DHI  •  DLOW  THEN  DIRFRAC  *  0  :  GOTO  8211 

8210  FOR  NREF  «  3  TO  87  STEP  3 

8211  IF  RNGREF(  NREF  )  <*  RMTR  GOTO  8220 

8212  RLOW  -  RNGREF(  NREF  -  3  )  :  RHI  -  RNGREF(  NREF  ) 

8213  IF  RLOW  -  RHI  THEN  REFFRAC  -  0  :  GOTO  8215 

8214  REFFRAC  *  (  RMTR  -  RLOW  )  /  (  RHI  -  RLOW  ) 

8215  IF  RLOW  >  0  GOTO  8217 

8216  RSHORT  -  BOUNDHI  «■  BOUNDLO  :  GOTO  8218 

8217  RSHORT  *  RAYREF(  NREF  -  3  );  RLONG  -  RAYREF(NREF) 

8219  REF  RAY  *  RSHORT  +  (  RLONG  -  RSHORT  )  *  REFFRAC  :GOTO 
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8400 

8220  IP  (  NREF  «  87  )  AND  (RNGREF(  87 )  >  0  )  GOTO  82  21  ELSE 
GOTO  8226 

8221  RLOtf  «  RNGREF(  87  ) 

8222  IF  RLOW  *  RMTR  THEN  REFRAY  =  RAY  REF  (  87  ):  GOTO  8400 

8223  REFFRAC  *  (  RMTR  -  RLOW  )  /  (  1 000  -  RLOW  ) 

8224  RSHORT  -  RAY REF (  NREF  )  :  ALONG  =  SQR  (  1000s 2  +  ( 
BOUNDHI  +  BOUNDLO  )s2  ) 

8225  REFRAY  «  RSHORT  +  (  ALONG  -  RSHORT  )  *  REFFRAC  :  GOTO 

8400 

8226  NEXT  NREF 

8227  GOTO  8410 

8370  REM  OUTPUT  FINAL  ANSWERS  FOR  EACH  DEPTH  AND  RANGE  TO 
DISK  FILE 

8385  LPRINT  "HYSPD-"  HYSPD 

8400  FINAL(  NMTR  )  =  HYSPD  /  (  REFRAY  -  DIRRAY  ) 

8410  NEXT  NMTR 

8420  MULT  »  ( TGTDPTH  l  10)  -  1 

8430  FOR  JMP  >1  TO  100 

8440  LSET  F1$  -  MKS$(FINAL( JMP)) 

8450  RECNUMt  ■  JMP  +  (MULT  *  1  00) 

8460  PUT  #1,  RECNUMft 

8470  NEXT  JMP 

8480  FOR  NUL  -  0  TO  100 

8490'  FINAL(NUL)  *  01 

8500  NEXT  NUL 

8510  FOR  NUL  »  0  TO  135 

8520  TEST(NUL)  *  01 

8530  NEXT  NUL 

8540  FOR  NUL  •  0  TO  15 

8550  NX( 15)  «  01 

8560  NEXT  NUL 

8570  FOR  NUL  »  0  TO  1 79 

8580RNGDIR(NUL)«01:RNGREF(NUL)=0:RAYDIR(NUL)=0: 

RAYREF(NUL) -0 
8620  NEXT  NUL 
9000  NEXT  NUPLMT 
9010  RETURN 

9500  PRINT:  PRINT:  PRINT  "INPUT  DATA  HAS  CAUSED  AN  ERROR" 

9510  GOTO  100 

9600  REM  PRINT  OUT  GRAPH  OF  HORIZONTAL  RANGE  .vs.  DEPTH 
9610  PRINT  CHR$(  26) 

9620  PRINT :PRINT:PRINT  "  SELECT  OUTPUT  DEVICE  PRINT 
9630  PRINT  1,"  VIDIO  SCREEN  ":  PRINT 
9640  PRINT  2,"  PRINTER  ":  PRINT 

9650  INPUT  "  ENTER  THE  NUMBER  FOR  OUTPUT  DEVICE  DESIRED 

" ,NUMB 

9660  IF  (  NUMB  <  1  OR  NUMB  >  2  )  THEN  PRINT  "IMPROPER  INPUT 
”:GOTO  9620 
9665  PRINT:  PRINT 

9670  PRINT  "  THE  OUTPUT  MATRIX  IS  PRESENTED  IN  10  METER 
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INCREMENTS  FOR  BOTH  " 

9680  PRINT  "  HORIZONTAL  RANGES  AND  FOR  DEPTHS  OF  10  METERS 
TO  200  METERS.  " 

9690  PRINT  "  HORIZONTAL  RANGES  FROM  THE  RECEIVER  TO  THE 
TARGET  ,  IN  " 

9700  PRINT  "  GROUPS  OF  100  METERS ,  INCREASE  FROM  LEFT  TO 
RIGHT,  " 

9710  PRINT  "  ie.  10  TO  100  METER  RANGES  WILL  BE  DISPLAYED 
FIRST,  FOLLOWED" 

9720  PRINT  "BY  110  TO  200  METERS  AND  SO  ON  UP  TO  1000 
METERS." 

9730  PRINT  "  DEPTHS  OF  0  TO  200  METERS  WILL  BE  DISPLAYED 
FOR  EACH  GROUP  OF  " 

9740  PRINT  "  RANGES,  WITH  THE  FIRST  ROW  REPRESENTING  10 
METERS,  THE  " 

9750  PRINT  "  NEXT  ROW  20  METERS,  ON  DOWN  THE  PAGE  UNTIL  ROW 
20  IS  " 

9760  PRINT  "  PRINTED  FOR  A  DEPTH  OF  200  METERS.": PRINT: PRINT 

9770  FOR  I  »  0  TO  9 

9780  THIS  «  10  +  100  *  I 

9790  THAT  »  1 00  +  1 00  *  I 

9800  IF  NUMB  «  1  THEN  PRINT: 

PRINT  "GROUP  OF  RANGES  FROM  "THIS"  TO 
”THAT"METERS" : GOT O  9  820 
9805  LPRINT: LPRINT: 

9810  LPRINT: LPR INT  "GROUP  OF  RANGES  FROM  "THIS"  TO  "THAT" 
METERS" 

9820  FOR  J  «  0  TO  19 
9830  DPH  »  10  +  10  *  J 

9840  IF  NUMB  »  1  THEN  PRINT:  PRINT  "DEPTH  OF  "DPH" 

METERS":  GOTO  9860 

9850  LPRINT: LPRINT  "DEPTH  OF  "DPH"  METERS" 

9860  FOR  K  «  1  TO  10 

9870  RECNUM%  «  K  ♦  100  *  J  +  10  *  I 

9880  GET  #1 , RECNUM % 

9890  AOUT  »  CVS(  FI  $  ) 

9900  IF  NUMB  «  2  GOTO  9930 
9910  PRINT  USING  "####.#  ";AOUT; 

9920  GOTO  9940 

9930  LPRINT  USING  "####.#  "jAOUTy 

9940  NEXT  K 

9950  NEXT  J 

9960  NEXT  I 

9965  LPRINT 

9970  RETURN 

10000  REM  CLOSE  ALL  FILES  AND  EXIT  PROGRAM 

10010  PRINT  CHR$(  26): PRINT  "DO  YOU  REALLY  WANT  TO  QUIT" 

10020  INPUT  "ENTER  A  Y  FOR  YES,  AN  N  FOR  NO  ",Z$ 

10030  IF  Z$  «  "Y"  THEN  END  ELSE  RETURN 
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